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Abstract 

One  of  the  remaining  challenges  involved  in  modern 
naval  ship  design  and  analysis  is  to  account  for  the  effects 
of  breaking  waves,  spray  and  air  entrainment  on  the 
performance  and  non- acoustical  signature  of  a surface 
ship.  The  near  field  flow  about  a surface  ship  is 
characterized  by  complex  physical  processes  such  as:  (i) 
spray  sheet  and  jet  formation;  (ii)  strong  free-surface 
turbulence  interactions  with  (large-amplitude)  breaking 
waves;  (Hi)  air  entrainment  and  bubble  generation;  and 
(iv)  post-breaking  turbulence  and  dissipation.  The 
challenges  associated  with  this  task  are  twofold.  The  first 
is  robustly  simulating  the  large-scale  problem  which 
involves  the  flow  about  an  entire  surface  ship.  The 
second  is  the  development  of  physics-based  closure 
models  for  steep  breaking  waves  in  the  presence  of 
turbulence.  To  wit,  a two-pronged  approach  consisting  of 
developing  an  understanding  for  closure  model 
development  and  applying  cutting-edge  computational 
capabilities  has  been  developed  to  accurately  simulate  the 
free-surface  flow  around  naval  combatants. 

Using  high-resolution  direct  numerical  simulation  of 
the  Navier-Stokes  equations  employeeing  the  level  set 
method,  we  have  successfidly  simulated  an  ensemble  of 
unsteady  breaking  waves  at  Reynolds  numbers  0(1 03  4). 
This  includes  steady  and  unsteady  as  well  as  spilling  and 
plunging  events.  This  dataset  is  continually  being 
improved  upon  in  terms  of  depth  and  breadth  as  a direct 
result  of  this  Challenge  Project.  The  goal  of  this  core 
research  area  is  to  develop  understanding  of  the  physics 
of  breaking  waves  to  help  guide  the  development  of 
physics-based  breaking  wave  modes.  The  dataset  is  being 


used  for  the  evaluation  of  closure  models  for  inclusion  in 
current  larger  scale  simulations  such  as  large  eddy 
simulation  and  Reynolds-Averaged  Navier-Stokes. 

Robustly  simulating  the  near-field  flow  of  a surface 
ship  requires  the  development  of  new  models  and 
numerical  techniques  suitable  for  use  in  large  scale 
applications.  We  have  performed  more  moderate-scale 
simulations  to  design,  verify,  and  validate  these 
capabilities  before  their  implementation  on  the  large- 
scale  simulations. 

Using  Numerical  Flow  Analysis  (NFA),  simulations 
of  several  naval  combatants  were  performed  at  a range  of 
speeds.  The  numerical  results  show  wave  overturning  at 
the  bow  and  flow  separation  at  the  transom.  Air  is 
entrained  along  the  side  of  the  hull  and  in  the  rooster-tail 
region  behind  the  stern.  In  both  regions,  numerical 
predictions  agree  well  with  experimental  measurements. 
This  work  marks  the  first  time  that  NFA  has  been  used  to 
simulate  an  entire  ship  hull.  The  numerical  simulations 
were  performed  on  the  Engineer  Research  and 
Development  Center  (ERDC)  Cray  XT3  using  128-256 
processors.  Approximately,  90  million  grid  points  were 
used  in  the  simulations. 

1.  Introduction 

The  Navy’s  design  requirements  for  faster  ships 
operating  in  unconventional  settings  are  pushing  the 
envelope  of  current  modeling  capabilities.  The  predictive 
uncertainties  introduced  by  breaking  waves,  spray 
formation  and  air  entrainment  have  ramifications  on  the 
signatures  of  proposed  vessels  that  in  turn  limit  the 
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confidence  in  analysis  of  new  designs.  The  Navy  is 
attempting  to  incorporate  high  performance  predictive 
capabilities  in  their  suite  of  design  tools  that  can  reliably 
predict  these  key  features.  This  Challenge  Project- 
supported  work  adopts  a two-pronged  approach  for 
developing  these  tools  by  advancing  the  state  of  the  art  in 
general  by  providing  robust  solvers  for  complex  fluid 
flows  while  concurrently  developing  physics-based 
models  for  complex  breaking  wave  and  air  entrainment 
phenomena. 

A simple  engineering  analysis  demonstrates  that 
these  flow  phenomena  depend  on  multiple  and  disparate 
time  and  length  scales.  At  the  smallest  scales,  spray 
formation  and  the  bubbly  contact  line  are  dependent  on 
length  scales  orders  of  magnitude  smaller  than  the  length 
of  the  ship  with  correspondingly  small  time  scales  as  well. 
At  the  largest  scales,  the  overall  structures  of  the 
shipwake  including  the  unsteady  rooster  tail  and  breaking 
waves  along  the  hull  are  a function  of  the  beam  of  the 
ship,  ship  speed,  etc.  Capturing  the  effects  of  the  very 
small-scale  physics  on  these  large-scale  features  in  a 
practical  design  environment  is  a primary  concern  in 
Computational  Ship  Hydrodynamics.  Even  with  the 
progress  that  has  been  made  in  computational 
methodologies  and  the  increase  in  computational 
capabilities  available  to  the  user,  a single  simulation 
which  hopes  to  resolve  these  phenomena  would  be  a 
monumental  undertaking  and  could  not  be  completed  in 
the  near  future. 

There  is  hope  that  these  physics  may  be  simulated 
directly  using  massively  parallel,  large-scale  applications 
with  models  to  bridge  the  gap  between  physics  and 
resolvable  scales.  This  requires  a robust  flow  solver 
capable  of  accurately  simulating  complex  flows  with 
minimal  user  guidance  which  ultimately  incorporates 
physics-based  models.  The  Numerical  Flow  Analysis 
(NFA)  code  provides  this  turnkey  capability  to  simulate 
the  near-field  flow  around  a surface  ship.  The  Level  Set 
Large  Eddy  Simulation  (LS-LES)  provides  the  framework 
for  model  development. 

The  goal  of  this  work  is  to  use  high  resolution 
simulations  of  the  complex  physical  process  involved  at 
their  most  basic  level,  to  develop  necessary  models,  and 
to  incorporate  them  into  large-scale  ship  simulations. 
These  results  must  then  be  compared  with  experiments  to 
determine  if  the  models  were  sufficient.  The  overall 
result  will  be  a general,  robust  numerical  capability  for 
use  in  the  design  of  surface  ships  in  a practical  design 
environment. 

2.  Problems  and  Methodologies 

Two  methodologies  specifically  designed  for  the 
scales  involved  are  used.  NFA  uses  Cartesian-grid  and 


volume-of-fluid  methods  to  simulate  the  flow  around 
surface  ships.  The  Cartesian-grid  portion  of  the 
formulation  ensures  that  complications  due  to  gridding 
are  eliminated.  The  sole  geometric  input  into  NFA  is  a 
surface  panelization  of  the  ship  hull  which  is  used  to 
impose  no-flux  boundary  conditions  on  the  hull  surface. 
The  volume-of-fluid  portion  of  the  numerical  algorithm  is 
used  to  capture  the  free-surface  interface,  including  large- 
scale  effects  of  breaking  waves,  spray  formation,  and  air 
entrainment.  The  Euler  equations  are  solved  using  a 
domain  decomposition  method.  Processor 

communication  on  the  Cray  XT3  is  performed  using 
either  Cray’s  shared  memory  access  library  or  Message 
Passing  Interface  (MPI).  The  CPU  requirements  are 
linearly  proportional  to  the  number  of  grid  points  and 
inversely  proportional  to  the  number  of  processors. 
Together,  the  ease  of  input  and  usage,  the  ability  to  model 
and  resolve  complex  free-surface  phenomena,  and  the 
speed  of  the  numerical  algorithm  provide  a robust 
capability  for  simulating  the  free-surface  disturbances 
near  a ship.  Additional  details  of  the  NFA  formulation 
are  provided  in  Dommermuth,  et  al.,  2006. 

A complementary  version  of  the  NFA  code  was 
developed  at  the  Massachusetts  Institute  of  Technology 
(MIT)  for  use  as  a test  bed  for  introducing  models  and 
improvements  to  NFA.  This  version  of  the  code-  is 
written  in  FORTRAN95  and  is  completely  vectorized  to 
maximize  performance  on  high  performance  computing 
(HPC)  platforms.  In  addition  to  these  coding 
improvements,  newly  developed  algorithms  are 
implemented  which  are  formulated  in  a general  way  to 
minimize  conditional  statements.  Examples  are  the  use  of 
a parametric  representation  for  the  ship’s  hull  and  general 
analytic  body  boundary  conditions.  These  improvements 
can  reduce  computational  costs  by  an  order  of  magnitude 
in  key  subroutines.  Further  information  regarding  these 
algorithms  can  be  found  in  Weymouth,  et  al.,  2006. 

The  small-scale  application  must  have  many  of  the 
features  of  the  large-scale  problem.  It  must  robustly 
handle  interface  topology  changes  without  ad-hoc 
treatments  or  human  intervention  and  have  very  little 
numerical  dissipation  which  may  mask  physical 
dissipation.  This  feature  is  highly  desirable  for 
developing  turbulence  models  of  the  micro-scale  physics. 
Until  physics-based  models  are  developed,  Direct 
Numerical  Simulation  (DNS)  of  the  Navier-Stokes 
equations  are  required.  Once  models  are  developed,  LES 
will  be  used  to  further  our  knowledge  of  larger-scale 
problems.  To  meet  this  goal,  LS-LES  was  developed  by 
the  Massachusetts  Institute  of  Technology.  It  is  a 
Cartesian  grid,  Eulerian-interface-capturing  scheme  based 
on  a modified  Level  Set  (LS)  method.  Both  water  and  air 
are  treated  as  incompressible  fluids  governed  by  a multi- 
fluid  Navier-Stokes  equation.  The  interface  between  the 
two  fluids  is  smoothed  over  a few  grid  points  for 
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numerical  stability.  The  Cartesian  grid  uses  a staggered 
MAC-type  formulation;  thus  no  upwinding  of  the 
convective  terms  is  necessary.  Processor  communication 
is  accomplished  via  MPI  making  the  code  portable  to  any 
platform.  In  particular,  this  code  performs  well  on  the 
IBM  P4,  the  Cray  XI,  and  the  Cray  XT3.  The  CPU 
requirements  are  linearly  proportional  to  the  number  of 
grid  points  and  inversely  proportional  to  the  number  of 
processors.  Both  large-scale  and  small-scale  applications 
have  been  validated  against  experiments  and  optimized 
for  use  on  HPC  resources  through  the  support  of  previous 
Office  of  Naval  Research  (ONR)  programs  and  Challenge 
Projects.  Additional  details  of  LS-LES  are  provided  in 
Hendrickson,  2004. 

3.  Results  & Discussion 

3.1.  LS-LES 

Flow  behind  a submerged  object  produces  a wave 
train,  one  to  two  body  lengths  behind  the  body.  Providing 
the  object  is  near  enough  to  the  surface,  the  first  wave  of 
the  resulting  wave  train  will  form  a quasi-steady  breaker 
(Duncan,  1983).  DNS  of  quasi-steady  breaking  waves  is 
performed  by  simulating  the  flow  past  a submerged  object 
of  diameter  1 centered  at  x=0.  An  immersed  boundary 
condition  method  is  used  to  simulate  the  presence  of  a 
submerged  cylinder  using  LS-LES.  The  key  focus  of  this 
work  is  to  use  the  geometry  to  generate  a quasi-steady 
breaking  wave.  Thus,  viscous  effects  about  the  body  are 
not  considered  and  the  condition  enforced  at  the  body  is  a 
no-penetration  condition.  Figure  1 shows  a waterfall  plot 
of  the  surface  elevation  for  the  latest  quasi-steady 
breaking  wave  simulation  for  the  wavebreaking  database. 
Included  also  are  contours  of  the  transverse  vorticity  at 
the  air-water  interface.  The  flow  is  from  the  left  to  the 
right  while  the  body  is  fixed.  Time  increases  along  the 
vertical  profile.  This  simulation  is  for  a moderate 
Reynolds  number  of  800.  Surface  tension  effects  are 
included  (We=puD2/a=10)  and  can  be  seen  in  the 
presence  of  capillary  waves  on  both  wave  crests,  evident 
by  the  bands  of  vorticity  shown  in  the  plot.  The  motion 
of  the  capillary  waves  on  the  front  of  the  first  wave  crest 
between  1 < t < 2 represent  an  initial  development  of  the 
wave  due  to  the  manner  in  which  the  body  is  slowly 
introduced  into  the  uniform  inflow.  The  primary  wave  of 
interest  is  the  second  crest  in  the  image,  first  behind  the 
body.  As  it  forms  a bulge  and  capillary  wave  train 
system,  it  becomes  unstable  and  the  wave  breaks  by  the 
bulge  moving  down  the  face  of  the  wave.  However, 
because  of  the  continual  source  of  energy  to  the  flow  from 
the  presence  of  the  body,  the  bulge  system  begins  to  set 
back  up  again  prior  to  the  end  of  this  simulation. 


A time  history  of  the  vorticity  field  in  both  the  air  and 
water  flows  during  this  event  is  shown  in  Figure  2.  Prior 
to  the  motion  of  the  bulge  (Figure  2a),  the  vorticity  field 
beneath  the  bulge  is  essentially  a single  (negative)  vortex 
beneath  the  bulge.  Near  the  onset  of  motion  of  the  bulge, 
the  presence  of  a small  positive  vortex  can  be  detected 
(Figure  2b,  x/D~3.2).  During  the  bulge  motion,  the 
vorticity  field  becomes  richer  in  that  small  vortex-pair 
structures  can  be  identified  as  well  as  the  presence  of 
large-scale  vortical  structures  in  the  wake  of  the  breaking 
crest.  Figure  3 shows  a close-up  of  the  bulge  region  after 
the  event  has  started.  The  cyclic  motion  of  the  bulge  on 
the  front  face  of  the  wave  has  been  seen  in  experiments 
(Duncan,  1983).  For  the  development  of  physics-based 
models,  a longer  simulation  is  required  to  catch  the  cyclic 
breaking  of  the  bulge.  However,  this  type  of  information 
is  critical  for  the  validation  and  development  of  physics- 
based  mixed-phase  turbulence  models  for  LES  and  “off- 
the-shelf’  Reynolds-averaged  Navier-Stokes  capabilities. 

3.2.  ATHENA  in  Waves 

Predictions  of  the  horizontal  and  vertical  loads  on  the 
model  hull  form  of  the  US  Navy  Research  Vessel 
ATHENA  held  fixed  in  regular  head  sea  waves  were 
made  with  the  NFA.  The  three-dimensional  numerical 
simulation  was  performed  using  850x192x128= 
20,889,600  grid  points,  5x8x4=160  sub-domains,  and  160 
nodes  on  the  Cray  XT3  at  the  ERDC.  The  length,  width, 
depth,  and  height  of  the  computational  domain  are 
respectively  4.0,  1.0,  1.0,  0.5  ship  lengths  (L0).  Grid 
stretching  is  employed  in  all  directions.  The  smallest  grid 
spacing  is  0.0020  L0  near  the  ship  and  mean  waterline, 
and  the  largest  grid  spacing  is  O.O20  L0  in  the  far  field. 
The  Froude  number  is  Fr  = Ua  / y]gL0  = 0.4136  where  U0 
is  the  ship’s  speed  and  g is  the  acceleration  of  gravity. 
Two  incident  wavelengths  (A)  are  considered:  2JL0  = 2 
and  AJLa  = 1/2.  In  both  cases,  the  wave  steepness  is  H//l= 
0.06,  where  FI  is  the  wave  height.  The  equations  for  the 
wavemaker  are  imposed  ahead  of  the  ship.  Initial 
transients  are  minimized  by  slowly  ramping  up  the  free- 
stream  current  and  the  incident  wave  amplitude.  For  this 
simulation,  the  non-dimensional  time  step  is  t=0.0005. 
The  numerical  simulation  runs  10,000  time  steps 
corresponding  to  5 ship  lengths.  Each  simulation  requires 
about  80  hours  of  wall-clock  time. 

Tests  of  the  model  hull  form  were  done  at  the  David 
Taylor  Model  Basin  (DTMB)  in  Carderock,  MD.  Figures 
4 and  5 show  perspective  views  of  the  predicted  free- 
surface  elevations  around  the  ship  at  time  t = 5.  The 
predicted  and  measured  free-surface  elevations  as  a 
function  of  time  at  x = 1.5  are  shown  in  Figures  6 and  7. 
The  predicted  free-surface  elevations  are  ramped  up  to 
their  full  height.  This  minimizes  transients  associated 
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with  starting  up  the  numerical  wavemaker.  The  measured 
free-surface  elevations  show  slight  irregularities  that  are 
associated  with  limitations  with  the  wavemaker  at  DTMB. 

Figures  8 and  9 show  the  predicted  drag  force 
compared  to  measurements.  In  the  case  of  the  numerical 
simulations,  the  drag  is  initially  zero  because  the  model  is 
ramped  up  to  full  speed  from  zero  forward  speed.  The 
primary  harmonics  that  are  evident  in  the  plots  are  due  to 
the  incident  wave  forces.  We  speculate  that  the  higher 
harmonics  that  are  evident  in  the  laboratory  results,  which 
are  especially  evident  for  the  longer  wave  case,  are  due  to 
vibrations  in  the  structure  that  is  used  to  restrain  the 
model.  In  general,  numerical  predictions  and  laboratory 
measurements  agree  well. 

Figures  10  and  11  show  the  predicted  vertical  force 
compared  to  measurements.  The  displacement  has  been 
subtracted  out  from  the  results.  As  the  model  ramps  up  to 
full  speed,  a mean  suction  force  is  induced  on  the  model. 
The  oscillatory  portion  of  the  force  is  dominated  by 
hydrostatics.  Once  again,  numerical  predictions  and 
laboratory  measurements  are  in  good  agreement. 

It  is  not  possible  to  show  all  of  the  details  in  the 
numerical  predictions  through  the  use  of  figures.  In  order 
to  study  the  flow  in  even  more  detail,  several  animations 
have  been  prepared  at  the  flow  visualization  center  at 
ERDC.  The  animations  show  forced  motions  of  naval 
combatants,  naval  combatants  in  waves,  and  naval 
combatants  moving  with  constant  forward  speed.  The 
animations  are  viewable  at  http://www.saic.com/nfa. 
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Figure  1.  Waterfall  plot  of  surface  elevation  (colored 
by  vorticity)  for  quasi-steady  breaking  wave.  The 
surface  profiles  are  increased  by  a distance  of  0.1  for 
each  T=1.0. 


Figure  2.  Evolution  of  vorticity  structure  in  both  air 
and  water  fluids  for  a quasi-steady  breaking  wave,  (a) 
prior  to  bulge  motion;  (b)  near  onset  of  bulge  motion; 
(c)  bulge  motion;  (d)  bulge  motion,  with  large-scale 
structures;  and  (e)  bulge  beginning  to  set  back  up. 
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